The Kolmogorov turbulence theory in the hght of 
six-dimensional Navier-Stokes' equation 



Shunichi Tsuge§ 

Abstract. The classical turbulence theory by Kolmogorov is reconsidered using 
Navier-Stokes' equation generalized to 6D physical-plus-eddy space. Strong pseudo- 
singularity is shown to reveal itself along the boundary 'ridge' line separating the 
dissipation and inertial sub-ranges surrounding the origin of the eddy space. A 
speculation is made that this singularity is generated by two dipoles of opposite sign 
aligned on the common axis. It is supported by the observation that the universal power 
spectrum calculated rediscovers the Kolmogorov's -5/3 power law as independent of 
the dimensional approach. 



1. Introduction 

As early as in 1941 Kolmogorov|l, predicted some universal features of fluid turbulence 
which have been confirmed by experiments in later years. It is rather surprising that 
they are derived from dimensional analysis based on the simple assumption of 'local 
homogeneity' for small-scale turbulence. More surprising is the fact that those have 
little to do with the equation of dynamics of fluid. 

It is addressed in this paper to locate this 'missing link' through rederivation of the 
-5/3 power law of the spectrum for the inertial subrange as the universal law using the 
dynamical equation proposed in ref.2. The key issue for this task to work is that equation 
on which to describe the features of small eddies as independent of individual flow 
geometry. To meet this purpose a six-dimensional Navier-Stokes equation is employed, 
where additional 3D space having a length dimension , corresponding to eddy size, is 
introduced. 

Originally such an equation has been derived elsewhere^, using non-equilibrium 
statistical mechanics starting from Liouville's equation. In this formalism the 3D 'eddy' 
space has been introduced as natural consequence of a mathematical procedure, namely, 
the separation of variables of turbulent fluctuation-correlation equation. The equation 
stands as 6D generalization of the Navier-Stokes equation which in 3D physical space 
degenerates to the classical equation. 

In what follows, however, it is intended to rederive the same equation using 
phenomenologies alone on the basis same as the classical theory of Karman and 
HowarthjS]. Then, the wave-number space is introduced for separating variables 
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(Sec. 2). It is then Fourier-transformed into eddy space, thereby the 6D Navier-Stokes 
equation, the basis to all what follows is established (Sec. 3). The equation gives a novel 
expression for turbulent dissipation, enabling to predict existence of a pseudo-singularity 
surrounding the dissipation region (Sec. 4). In inertial subrange, this singularity turns 
out to be a couple of dipoles of opposite sense, separated by the order of several tens to 
hundreds the Kolmogorov length (Sec. 5). Power spectrum of this locally homogeneous 
turbulence is calculated, showing a wave- number dependence close to -5/3 power law of 
the Kolmogorov theory. (Sec. 6) 



2. Karman-Howarth formalism revisited as governing inhomogeneous 
turbulence 

In 1938, von Karman and Howarth proposed an equation governing homogeneous 
isotropic turbulence whose original form is written asf^ 



m;.(ns), + mKns), = (1) 

In this equation u'j = u'j{x, t) and u'l = u[{x, t) are instantaneous velocity fluctuation at 
X and X , respectively, overbar ( ) denotes the conventional (ensemble) average and 

NS = NS(V,n,p) 

d 

= (^+M- V-z/V^)w + p~Vp = (2) 

denotes the Navier-Stokes equation written in terms of the instantaneous fluid quantities 

u = u + u' 

p = p + p 

where u and p are the average velocity and pressure, respectively, u is the kinematic 
viscosity, and V denotes nabla vector. Likewise NS is defined by Eq.(2) in which ( 
V,u,p ) is replaced with ( V,n, p). 

Eq.(l) is an equation in 6-D physical space and time, namely, has seven independent 
variables {x,x,t). The classical theory jH] chose to discuss homogeneous and isotropic 
turbulence to reduce independent variables to {\x — x\,t). We can show, however, that 
the equation can be renovated to be available for inhomogeneous turbulence as well 
by introducing separation of variables, together with a closure condition which seems 
plausible intuitively. 

The actual derivation proceeds as follows: Eq.(l) in which NS(V,ii, p) is 
decomposed into average and fluctuating parts using (3), reads 



^(NS)l + wKNS); = (4) 

where 



together with (NS); defined as in Eq.(l). Upon substituting Eq.(5) into Eq.(4) we see 
that the equation consists of terms of double and triple fluctuation correlations, which 
we will decompose as 

gj{x, k)gx{x, k)S{k + k)dkdk , X^{1, 4) (6) 

-oo 

/oo „ _ „ _ „ _ 

gj{x, k)gi{x, k)gr{x, k)6{k + k + k)dkdkdk (7) 
-oo 

where k is the wave number, / is the characteristic length of flow geometry, R.P. denotes 
taking real part, and subscript A runs 1 through 4. A supplementary definition 

< = P~ V (8) 

stands for the pressure fluctuation. It is more convenient for later use to employ the 
following alternative expressions: 

/oo 
gj{x,k)gxix,-k)dk, A=(/,4) (6') 
-OO 

/OO roo ^ ^ ^ 

gj{x,k)dk / gi{x, —k + k)gr{x, —k)dk , {x;xorx) (7') 
-oo J —oo 

It is to be noted that there is no a priori reason why the double and triple 
correlations are to be expressed in terms of the same g\. At this point it is simply 
invoked as the closure condition to truncate the chain of equations at the level of Eq.(4). 
This closure can be justified only through the quality of outcome as to whether it fits 
with physical reality of turbulence. 

With these preliminaries, Eq.(4) with assumption (6 ',7') substituted into, leads to 
the form typical of separation of variables, here into (a;,t) and {x ,t), respectively; 

dHsilsr'Wf -i-9r'(rs)<°'] = o (9) 



/ 

J — ( 



= iuj = -^00 

where uj is the separation parameter having the dimension of the frequency, and gi 
stands for gi{x,k) , and (ns)*^°'' is defined by 



with 



/oo ^ ^ ^ 

gj{k - k)gr{k)dk 
-oo 

Imaginary factor i as appearing in Eq.(9) reflects the statistical symmetry of the tensor 
u'ju'i , which is met under the following additional condition 

gi{x,-k) = [gi{x,k)Y (11) 

In fact, then, the two terms inside[ ] of integral (9) are commutable to each other 
through taking complex conjugate (*), if iuj is purely imaginary. Thus we are led to the 
assertion that the only equation that need to be solved is 

lujgj = (ns)f (12) 
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To go further from this point on, we need to have relationship between frequency 
uj introduced as the separation parameter and wave number k by which the solution is 
constructed in the form of bilinear integral. Physically they are related to each other 
by the dispersion relation uj{k) . Or, instead, we may introduce phase velocity c by 

UJ = c - k (13) 

with no loss of generality. In fact, its constancy does not mean the turbulent eddies 
being assumed as nondispersive. For the unsteady term in the equation {d/dt ^ 0) is 
responsible for the dispersive part, if any. 

Owing to the convolution form of integral (10) periodic factor drops off from Eq. (12) 
despite its nonlinear structure by putting 

g^{x,k) = e^-"" h{x,k) , \ = {IA) (14) 

thereby we have the following equation governing its amplitude /a, 

iujfj = (ns)j (15) 



where 



(ns)j = [ns(/,/4 



f) fyoi ■ 

[_ + urdrik) - Ud^mfj + + di{k)U + dr{k)I{f,fr) (16) 



with 



d{k)=Vix)+ik (17) 

It is readily checked that expressions defined by (15) and (10) are related to each other 
by 

(ns)? = [(ns),]a(fc)^v(x) 
Several remarks are in order with regards to physical implications of variables 
having appeared in this section: Amplitude function / as obeying Eq.(14) is essentially 
complex and is not an observable. It is related to the observable quantity of fluctuation 
correlation through 

^ = R.P./=^ r e'^-(^-*Vj(a;, k)f;{x, k)dk (18) 

J —oo 

as is confirmed by (6'),(ll)and (13). Its complex variable structure has its origin 
in Eq.(9) where the imaginary unit i is introduced from the symmetry postulate of 
statistical mechanics j2] for correlation tensor (18). We may note some coincidental 
parallelism to Schroedinger's wave equation where the imaginary factor secures the 
corresponding tensor to be Hermitean as it should. 

The most crucial on which this paper rests is the soundness of the basis of Eq.(l). 
Originally it was a product of intuition by Theodore von Karman without any 'first 
principle' ground shown in the paper In fact, then, any linear combinations of 
(NS)' would be claimed as equally qualified. Firm basis is provided by nonequilibrium 
statistical mechanics^!] to justify that the linear combination only of the form (1) is 
consistent with Liouville's equation, namely, the equation of continuity in the phase 
space of Hamiltonian mechanics. . 
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3. Turbulence in eddy space 

The nonlinear integro-differential equation we have derived in the previous section is not 
very easy to deal with. Since, however, the integral is of convolution type, it transforms 
to a simple product, through Fourier transform as 



fj{x,k) = {27Tiy 



(19) 



The actual operation of the transform on Eq.(14) gives 

d d dii ■ 

[-Cr-^ + -^ + Urdr{s)-ud'^{s)]qj + dj{s)q4 + -^qr + dr{s){qjqr) = 0(20) 

where 

dj{s)^d/dx,+d/dsj (21) 

Variable s introduced in Eq.(19) as the Fourier variable adjoint to wave number k has 
dimension of length, which stands for eddy size with direction of the vorticity vector, 
so may well be called eddy variable. It should be remarked that Eq.(20) is the Navier- 
Stokes equation in 6D (physical plus eddy) space, describing the motion of turbulent 
vortices moving with phase velocity c . In fact, this equation governing qx{x,s) for 
prescribed u{x) and p{x) is alternatively written as 

[c • V(s)]g = NS(a(s), u + q,p + pq^) - NS(V(a;), w,p) (22) 

with 

d{s) = V(x) + V(s) (23) 

where NS has been defined by Eq.(2), and V denote the nabla vector in the respective 
spaces§. 

Simplicity in expression for physical quantities is another advantage of working with 
s-space. For example, fluctuation-correlation formula (18) takes remarkably simple form 

/oo 
dsqj(x,s + x)qi(x,s + x) (24) 
-oo 

as is easily confirmed by substituting (19) and employing definition of the delta function 
S{s/l) = {2niy J e'^'^dk . Turbulent dissipation which has been dealt with as an 
elementary parameter in the classical dimensional analysis can be expressed in this 
space by an explicit form: We have, by definition. 



e V f duj ^ du'i 
2 2 1 dxi dxi 



du'i _^ du', 



du'i 



1 du'i 

^ + — - 



dxi dxj I \ dxi dx 



(25) 



J / J x=x 



— V 



' &^u!jU'^^ 
dxidxi 



+ 



d^u'jU'i 
dxidxj 



J x=x 



The bold-face letters denote vectors having three components, not to be confused with 6D vectors. 
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On the other hand, we have from (24) 



dxmdxn {27rl)^ J-oo \dxm. ds 
1 

ds[dm{s)qj][dn{s)qi] 




s=s 



Thus we are led to the final expression for dissipation e as 



2 = j_Jsms)q,f + {di{s)q,){d,{s)qi)] (26) 

This formula will reveal a new facet of the dissipation function hidden in the classical 
working space. 

4. Dimensional consideration on the existence of pseudo-singularity 

In the Kolmogorov regime, dissipation e was introduced, together with kinematic 
viscosity v , characteristic velocity v and length / , respectively, as an elementary 
parameter by which to construct the dimensional analysis. They are related to each 
other through the following relationships 

e ~ v'^/l ~ 4//k (27) 

where the subscript K denotes the respective quantities in the dissipation region. The 
following formulae are their immediate consequences 

vkIk = y 

= vR-^l^ \ (28) 

Ik = IR~^^^ (Kolmogorov length) 

where R is the fiow Reynolds number 

R = vl/v (29) 

Formula (26) we have derived in the newly defined hyperspace sheds some lights on 
this classical theory that is conducted within the physical space, in the sense that e is not 
necessarily an elementary parameter any longer. According to Kolmogorov[T] viscous 
dissipation occurs only within the scale of Kolmogorov length. Then, the integral region 
of dissipation function (26) which is proportional to the kinematic viscosity is confined 
within a small volume of 0(/^)- Dimensional analysis of (26) gives 

'-^F©^'- (30) 

Thus, from (27) through (30) we are given for the order of magnitude of |q| : 

\q\lv ~ R"^ (31) 

It is an indicative of strong in/out mass fiow existent (cf. the second of Eq.(28) of the 
classical theory) in the vicinity of this extremely small 'energy black hole' where the fiow 
loses kinetic energy e converted into heat. Inside this region |s| < 0(/k) the turbulence 
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dies off rapidly towards g(0) = .Tliere must be, tlierefore, a drastic variation in tlie 
magnitude of |q| wliicli peaks at tlie boundary ridgeline between dissipation and inertial 
ranges surrounding tlie origin (s=0), and diminishes quickly inside. 

In the neighborhood of this 'pseudo'-singularity it is obvious for the following 
conditions to hold: 





(32) 

The latter condition warrants for local homogeneity in the sense of Kolmogorov to hold 
most strictly in the localized eddy space for any inhomogeneous turbulence. Under these 
circumstances, the 6D Navier-Stokes equation (20) together with equation of continuity 
reduce to 

(33) 

*^ - "5) * ^ i = " 

which is nothing but the classical (3D) equations for laminar viscous flows, as disguised 
through 

s , {u,p) {q,pq4) (35) 

This rule reigns dissipation region as well as the adjacent region of inertial subrange 
surrounding it, where q-function diminishes outward off the pseudo-singularity until 
flow inhomogeneity starts to make its appearance. The solution as such, to be pursued 
in this space, is 'universal', namely, to be valid even for any shear turbulence. 

According to Kolmogorov jl], the inertial range is where no viscous effects are 
operating. So we may claim that the potential flow is prevailing there. This speculation 
is supported by the assertion: ^ potential flow is the solution of the Navier-Stokes equation 
in the region where no solid boundary is existent.^ Now is the case with it, because no 
physical substances are intervening here in this space. 



5. The pseudo-singularity in the inertial subrange 

The pseudo-singularity, as viewed from the domain of the inertial range, looks as if a 
genuine singularity, whose actual form is now to be identified. 

We start with checking if the local isotropy in the classical sense has physical reality. 
The isotropy assumption is equivalent to 

qj = SjQ{s) (36) 

which is the alternative expression of Robertson's theoremj^.^- Substitution of (36) into 
equation of continuity (33) gives 

sdQ/ds + 3Q = 
.-. Q ~ 
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Obviously, this is a source/ sink flow, corresponding to the velocity potential 

= as-^ (37) 

If a is positive it represents a sink flow with mass flux Ana vanishing at the origin. 
Reminding that the dissipation range |s| < /k is the 'black hole' of turbulent kinetic 
energy losing the amount by e converted into heat every second, we see that no mass is 
supposed to be lost. Thus solution (36) for spherical isotropy is ruled out. 

Prospective singularity now to take over can be sought within potential flow regime 
as follows: It is obvious that operator / ds\ds^ds'^{N = I + m + n) commutes with 
Laplacean operator, therefore 

_ 9"^'" (38, 
represent a group of potential flows. In particular 



^ dsi 



(39) 

representing a dipole aligned with its axis parallel to Sj direction is seen to meet the 
purpose. In fact a dipole is made up of a pair of sink/source of equal strength, so no 
mass flux is lost at this spot. Thus simplest possible candidate for the pseudo-singularity 
is axially isotropic. 

The following observation from direct numerical simulation [H] helps us draw a more 
precise picture of our pseudo-singularity: In the physical space the dissipation occurs 
only within a conflned volume of 'elementary particles' of rod shape, randomly dispersed 
in turbulent medium. They all have a shape like a worm with radius (Tq of several 
Kolmogorov lengths and several ten times of it lengthwise ( sq = Nao ), rotating round 
their own axes. Their lifetime is about Ik/vk = {l/v)R~^^'^ . So to say, they are like flrefly 
worms illuminating light for their lifetime as small as of the order of far subseconds. 

This picture when mapped onto our eddy space is such that two spinning dipoles 
having opposite sign and rotation are placed at 

S= (S1,S2,S3) = (±So,0,0) (40) 

respectively, where the axis of symmetry is Si axis. The direction j = 1 is the direction 
of motion of turbulence-generating body, representing the only vector prescribing the 
fluid motion. The flow fleld induced by the pair of spinning dipoles 

q = Qd + Qs 

Qd = Qd + Qd \ (41) 
9s = Qs+Qs , 

with and standing for velocities induced by dipoles and line vortices, placed at 
points (40), respectively: 



9 



where i is the unit vector designating Si axis. Expression (42) is direct consequence of 
(37) and (39), and that for (43) is Biot-Savart's law for a hne vortex with infinitesimal 
length directing Si-axis and with circulation f3. For example, the actual form of qi{s) 
rewritten in axially isotropic form qi{si, a) with a"^ = $2 + is, 

d"^ ( a a ] 



dsl \ [(si - So)2 + a2]i/2 [{s, + Soy + (72]i/2 

3a2 



+ 



[(Si - So)2 + (72]3/2 [(s, _ So)2 + ^2]5/2 

2 3(72 ^ 



(44) 



[(Si + So)' + [(si + So)' + <J']'/' j 

Note that there is no contribution from qg to gi. For >> so this expression 
approaches to 

qi ~ 6asoSi H (45) 



(46) 



(47) 



which is a quadrupole field, as it should be expected. 

Streamlines are shown in Fig.l of the fictitious fiow generated by a pair of dipoles 
given by the potential 

0(1) ^ d(P^^^/dsi 

0(0) = oc{[{s,-Sof + a^]-^"-[{sr + Sof + a^Yy^} 
This axi-symmetric fiow can also be represented using streamfuntion '0 as 

= = a-^d{ai})/da 1 

= d(t)^^yda ^ -a-^d{aij)/ds^] 
from which we have 

^ = -(790(°)/9(7 (48) 

The fiow pattern :const. shows quadrupole- like structure at far field ( \s\ » Sq 
) toward which the longitudinal vortices are stretched streamwise and then getting 
thicker. On the returning path to the dipole core they are chopped ofi^ and trim the 
aspect ratio, getting into the dissipation region. This picture may serve the qualitative 
description of what is actually observed (Fig. 2). 

Expression (45) gives us estimate for the outer boundary of the locally homogeneous 
region. That is also the inner boundary of inhomogeneous region where the first 
condition of inequality (32) ceases to hold; 

|g| ~ V (49) 

Let this boundary be defined by s ~ /o, then we have q ~ asol^'^ from (45), therefore 
condition (49) is replaced with a more precise one 

vo ~ aso^o"*^ (50) 
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where Vq is the characteristic velocity corresponding to Iq. They supplement Kolmogorov 
formula (27) as 

e ~ v^/Ik ~ vl/lo ~ vyi (51) 

On the opposite side of the inertial subrange s — sq ~ /k, cr ~ ^k, 9 is estimated from 
(44) as 

q{so) ~ «/k', (52) 

We note that the dimensional analysis developed in Sec. 4 still holds by replacing \q\ 
with q{so) of (52), for instance, 

q{so)/v ~ R'/'. (53) 

Then, by eliminating a, q{so) and vq from (50) through (53) we have 

lo/l = {so/ir^''R-''^'' (54) 

The relationship between sq and Ik or / is yet to be reconsidered. At present no 
consensus formula is available for explicit parameter dependence of worm size sq ■ A 
recent observation by numerical experiments [H] is that sq is of the order of the Taylor 
microscale (~ lu'/v, u' ; r.m.s. of the velocity fluctuation ), according to which the 
outer boundary of the inertial range is 

lo/l ~ i?-ii/26(^/^)4/i3 (55) 

6. Power spectrum for inertial subrange 

The actual form of the pseudo-singularity as predicted in the preceding section enables 
us to calculate ID power spectrum for ID wave number ki in the inertial subrange. This 
spectrum function Pii(fci) is written, by definition. 



<= / Piiiki)dk, (56) 

J— oo 

The well-known consequence of the classical dimensional analysis is that Pu , having 
the dimension of v^l ~ e^/^/^^^, reads in the language of wave number space as 

Pn ~ fci"'^' (57) 

An alternative look into this formula on our own basis is the following: The actual 
form of Pu is obtained through comparing formula (18) for x = x and j = I = 1, 
namely, 

/oo /"GO 
dki / dk2dk3fif* 
-oo J —QO 

with (56), which reads 

/oo 
dk2dksfj* (58) 
-oo 

This integral is, upon substitution of (19), transformed into the one in s-space as 

/oo /"OO /"OO 

dsi dsie~'^^'^''''n ad(Tgi(si,a)gi(si,(T) (59) 
-oo J —oo Jo 
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Since the integration spans over the whole s-space, we have yet to know solution qi 
inside the dissipation range where essentially viscous flow prevails. However, 

the volume of the dissipation range is by far the smaller than inertial subrange [0(/o)], 
we may dispense with potential flow solution (44) to be integrated over its own region. 
Thus integral (59) with a small spheroidal regions excepted is to be carried out. (See 
Fig. 1.) Then, the integral is shown to be converted into a double integral as follows 

^^^^^^^^(2^/0 '^('^')^('^'' ^1)' (60) 
with Q(cr^, /ci) defined by 

roo ^20(0) 

Q{a ,ki)^ J sinkiSi—^-^-dsi (61) 
where boundary contour Si = Si(cr^) is given by. 
s\{a') ^ soil + S) 

The ID power spectrum calculated is shown in Fig. 3. Parameter 6 corresponds to 
the slenderness ratio of the worm as detected by the direct numerical simulation, which 
is estimated as the order of O(10iV)~^ , with N a number of the order of unity to several. 
For a certain range of this parameter, the spectrum shows ki dependence, then with 
increase in 5 it transits to /cjf^ for 5 >> 1 , where the pair of dipoles are regarded as a 
quadrupole asymptotically. 

7. Conclusions 

6D Navier-Stokes' equation [Eq.(22)] governing turbulence is specialized to the locally 
homogeneous dissipation/inertial ranges in the eddy space. This equation [Eq.(34)] has 
revealed existence of a (pseudo-)singular solution of shape like a volcano with a ridgeline 
separating the inertial subrange from the dissipation range inside which is a energy 
black hole of turbulence. This finding of anomalous hike of the velocity wave function 
((31) or (53)) in the eddy space has enabled us to rediscover the universal form of the 
power spectrum from the equation of fluid dynamics as distinct from the dimensional 
analysis of the classical theory. In fact, an analytical solution of one-parameter family 
defining the ridgeline contour includes a case that is close to the consensus wave number 
dependence of -5/3 power law by Kolmogorov. It is yet left to be answered to eliminate 
the parameter dependence. It will be achieved by replacing the inviscid solution valid 
only for inertial subrange with a prospective viscous solution that is uniformly valid 
throughout dissipation and inertial regions. 
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Figure Captions 

Figure 1. Pseudo-singularity formed by a pair of dipoles in the inertial subrange of 
the eddy space. The shaded part represents the dissipation region. [See Eq.(63).] 
Figure 2. Evolution of turbulent vortices visualized by a gigantic turbulent jet : Volcano 
eruption of Mt. Saint Helens [7]. 

Figure 3. ID power spectrum Pu plotted against ID wave number ki with slenderness 
ratio 6 as the parameter. The calculated spectrum is averaged over five neighboring 
points to make the comparison with -5/3 or -2 power laws easier. 
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